######
### Load libraries and set seed
######

set.seed(1492)

library(foreign)
library(arm) #this is leading to some sort of conflict
library(blme)
library(texreg)
library(Hmisc)
library(lme4)
library(RColorBrewer)
library(boot) 
library(stargazer)
require(mgcv)
require(ggplot2)
library(lmtest)

##################
## FUNCTIONS 
##################


vcovCluster <- function(
  model,
  cluster
)
{
  require(sandwich)
  require(lmtest)
  if(nrow(model.matrix(model))!=length(cluster)){
    stop("check your data: cluster variable has different N than model")
  }
  M <- length(unique(cluster))
  N <- length(cluster)           
  K <- model$rank   
  if(M<50){
    warning("Fewer than 50 clusters, variances may be unreliable (could try block bootstrap instead).")
  }
  dfc <- (M/(M - 1)) * ((N - 1)/(N - K))
  uj  <- apply(estfun(model), 2, function(x) tapply(x, cluster, sum));
  rcse.cov <- dfc * sandwich(model, meat = crossprod(uj)/N)
  return(rcse.cov)
}


#######################################
# Multiple plot function
#
# ggplot objects can be passed in ..., or to plotlist (as a list of ggplot objects)
# - cols:   Number of columns in layout
# - layout: A matrix specifying the layout. If present, 'cols' is ignored.
#
# If the layout is something like matrix(c(1,2,3,3), nrow=2, byrow=TRUE),
# then plot 1 will go in the upper left, 2 will go in the upper right, and
# 3 will go all the way across the bottom.
#
##########################################
multiplot <- function(..., plotlist=NULL, file, cols=1, layout=NULL) {
  require(grid)

  # Make a list from the ... arguments and plotlist
  plots <- c(list(...), plotlist)

  numPlots = length(plots)

  # If layout is NULL, then use 'cols' to determine layout
  if (is.null(layout)) {
    # Make the panel
    # ncol: Number of columns of plots
    # nrow: Number of rows needed, calculated from # of cols
    layout <- matrix(seq(1, cols * ceiling(numPlots/cols)),
                    ncol = cols, nrow = ceiling(numPlots/cols))
  }

 if (numPlots==1) {
    print(plots[[1]])

  } else {
    # Set up the page
    grid.newpage()
    pushViewport(viewport(layout = grid.layout(nrow(layout), ncol(layout))))

    # Make each plot, in the correct location
    for (i in 1:numPlots) {
      # Get the i,j matrix positions of the regions that contain this subplot
      matchidx <- as.data.frame(which(layout == i, arr.ind = TRUE))

      print(plots[[i]], vp = viewport(layout.pos.row = matchidx$row,
                                      layout.pos.col = matchidx$col))
    }
  }
}


#############
### Import and Process Data
#############


### CD Presidential results
library(dplyr)
cd_pres_results<-read.csv(file="cd_pres_results.csv")
cd_pres_results2<-cd_pres_results
cd_pres_results2$year<-2006
cd_pres_results3<-cd_pres_results
cd_pres_results3$year<-2008
cd_pres_results4<-cd_pres_results
cd_pres_results4$year<-2010
cd_pres_results5<-cd_pres_results
cd_pres_results5$year<-2012
cd_pres_results<-rbind(cd_pres_results2, cd_pres_results3)
cd_pres_results<-rbind(cd_pres_results, cd_pres_results4)
cd_pres_results<-rbind(cd_pres_results, cd_pres_results5)
cd_pres_results$pres_dem<-NA
cd_pres_results$pres_dem[cd_pres_results$year==2006]<-cd_pres_results$ssp_dem_2004[cd_pres_results$year==2006]
cd_pres_results$pres_dem[cd_pres_results$year==2008]<-cd_pres_results$ssp_dem_2004[cd_pres_results$year==2008]
cd_pres_results$pres_dem[cd_pres_results$year==2010]<-cd_pres_results$ssp_dem_2008[cd_pres_results$year==2010]
cd_pres_results<-select(cd_pres_results, cd_fips, year, pres_dem)

load(file="spatial_voting_data.RData")

dwdime_rf<-read.dta(file="dwdime_rf.dta")
merged<-merge(merged, dwdime_rf, by.x=c("fips", "year"), by.y=c("fips", "cycle"), all.x=T)
### Analysis using Adams et al data

### Setup what kind of analysis we want to do.  This switch could be any of the four "analysis" types below.  
analysis<-"house_vote"
today<-format(Sys.time(), "%Y%m%d")

library(foreign)

data <- read.dta("adamsetal_operational ideology.dta")

summary(data$d_cand_ip)
summary(data$d_cand_ip[data$incumbent_party == -1])

summary(data$r_cand_ip)
summary(data$r_cand_ip[data$incumbent_party == 1])

normal_range <- data$d_cand_ip < 0.2 & data$r_cand_ip > 0.3 &
                data$r_cand_ip < 1.62

data$cutpoint<-(data$r_cand_ip-data$d_cand_ip)/2

data$cutpoint<-(data$r_cand_ip-data$d_cand_ip)/2
data$district<-paste(data$st, data$cd, sep="-")

#data$survey_ideology<-data$np_score
#data$std_survey_ideology<-scale(data$survey_ideology)
### Main Model
mod2 <- lmer(rhousevote ~ scale(np_score)*scale(cutpoint) + factor(pid7) + (1|district),data=data[normal_range,])
summary(mod2)



########################################
### Start analysis
#######################################


########################        
## The remaining code generates the main ggplot graphs.
######################## 

## This codes adds some variables we'll need to generate the graphs
merged$weight[is.na(merged$weight)]<-1
merged$weight<-as.numeric(as.vector(merged$weight))
merged$bin<-NA
## Break people into quintiles based on their ideology
ideo_quantiles<-quantile(merged$survey_ideology, na.rm=T,probs = seq(0, 1, 0.33))

## use weighted quantiles
ideology<-merged$survey_ideology[!is.na(merged$survey_ideology)]
ideo_quantiles2<-Hmisc::wtd.quantile(merged$survey_ideology, weights=merged$weight, probs=seq(0, 1, 0.33))

for (i in 1:3){
	merged$bin[merged$survey_ideology < ideo_quantiles2[i+1] & merged$survey_ideology > ideo_quantiles2[i]]<-i
	print(i)
}

merged$bin2<-NA
merged$bin2[merged$bin==1]<-"#3300FF"
merged$bin2[merged$bin==2]<-"#3333FF"
merged$bin2[merged$bin==3]<-"#3366FF"
merged$bin2[merged$bin==4]<-"#3399FF"
merged$bin2[merged$bin==5]<-"#33FFFF"
merged$bin2[merged$bin==6]<-"#FF6600"
merged$bin2[merged$bin==7]<-"#FF3300"
merged$bin2[merged$bin==8]<-"#FF0000"
merged$bin2[merged$bin==9]<-"#CC0000"
merged$bin2[merged$bin==10]<-"#990000"
 
colors2<-c("1"="#0003FF",  "2"="#3366FF","3"="#000000", "4"="#FF3300", "5"="#CC0000") 

merged$congress_post_vote[merged$congress_post_vote>1]<-NA
merged$congress_post_vote<-as.numeric(as.vector(merged$congress_post_vote))

merged$leg_pid3_2<-as.factor(merged$leg_pid3)
merged$pid3_2<-as.factor(merged$pid3)

merged<-merged[!is.na(merged$pid3),]
merged<-merged[!is.na(merged$leg_pid3),]
merged$leg_pid3<-factor(merged$leg_pid3)
merged$pid3<-factor(merged$pid3)

merged<-merged[merged$leg_pid3!=2,]

merged$pid3<-as.vector(merged$pid3)
merged$pid3<-factor(merged$pid3)

merged$leg_pid3_names<-as.vector(merged$leg_pid3)
merged$leg_pid3_names[merged$leg_pid3==1]<-"Legislator-Democrat"
merged$leg_pid3_names[merged$leg_pid3==3]<-"Legislator-Republican"
merged$leg_pid3_names<-factor(merged$leg_pid3_names)


line.type<-c( "longdash","dashed", "dotdash", "dotdash", "dotted")

ideo.col<-c("1"="#0003FF",  "2"="#3366FF","3"="#000000", "4"="#FF3300", "5"="#CC0000") 


############################################################
## This loop outputs the full panel of 6 loess curves; simpler graphs are below 
## (Figure 2 in Spatial Voting paper, October 2014)
############################################################

## Note for CW: I need to add comments in this chunk of code.

graph_types<-c("all_contested_races")
i<-1
merged2<-merged


x.text<-"Pr(Vote for Incumbent)"
range.dems<-c(-.6,-.1)  ##range of 90% of legislators, 5th percentile - 95th percentile
breaks.dems<-c(-.6,-.5,-.4,-.3,-.2,-.1)  ##range of 90% of legislators, 5th percentile - 95th percentile
range.reps<-c(.35,.9)  ##range of 90% of legislators, 5th percentile - 95th percentile
breaks.reps<-c(.4,.5,.6,.7,.8,.9)  ##range of 90% of legislators, 5th percentile - 95th percentile

#range.dems<-c(1,4)  ##range of 90% of legislators, 5th percentile - 95th percentile
#breaks.dems<-c(1,2,3,4)  ##range of 90% of legislators, 5th percentile - 95th percentile
#range.reps<-c(4,7)  ##range of 90% of legislators, 5th percentile - 95th percentile
#breaks.reps<-c(4,5,6,7)  ##range of 90% of legislators, 5th percentile - 95th percentile

merged2$dv<-merged2$depvar
merged2$iv<-merged2$leg_ideology

merged2<-merged2[merged2$incumb %in% c(-1,1) & !is.na(merged2$challengers.total.receipts),]

pid3<-c(1,2,3)

allvotes_demleg <- ggplot(merged2[ merged2$leg_pid3==1 & !is.na(merged2$weight),],  aes( x = iv, y = dv))
allvotes_demleg <- allvotes_demleg+ylab(x.text)+ xlab("Incumbent Position") + theme(legend.position="none")+scale_x_continuous(limits=range.dems, breaks=breaks.dems)+scale_y_continuous( breaks=c(0,.2,.4,.6,.8,1))+ ggtitle("All Voters/Democratic Legislator")  +coord_cartesian(ylim=c(0,1))+ theme(plot.title = element_text(size = rel(1))) + theme (axis.title.y  = element_text(size= rel(.75)))+ theme(plot.title = element_text(hjust = 0.5))


for (p in pid3){
	allvotes_demleg <- allvotes_demleg + stat_smooth(se=TRUE, method="gam", formula = y ~ s(x, k=4), aes(weight=weight), data=merged2[merged2$leg_pid3==1 & merged2$pid3==p & !is.na(merged2$weight),], colour="black",linetype=line.type[p], size=1)
}

allvotes_repleg <- ggplot(merged2[ merged2$leg_pid3==3 & !is.na(merged2$weight),],  aes( x = iv, y = dv))
allvotes_repleg<-allvotes_repleg+ylab(x.text)+ xlab("Incumbent Position") + theme(legend.position="none")+ scale_x_continuous(limits=range.reps, breaks=breaks.reps)+scale_y_continuous( breaks=c(0,.2,.4,.6,.8,1))+ ggtitle("All Voters/Republican Legislator")  +coord_cartesian(ylim=c(0,1)) + theme (axis.title.y  = element_text(size= rel(.75)))

for (p in pid3){
	allvotes_repleg <- allvotes_repleg + stat_smooth(se=TRUE, method="gam", formula = y ~ s(x, k=4), aes(weight=weight), data=merged2[ merged2$leg_pid3==3 & merged2$pid3==p & !is.na(merged2$weight),], colour="black",linetype=line.type[p], size=1,linetype=line.type[p], size=1)+ theme(plot.title = element_text(size = rel(1),hjust = 0.5))

}

leg_pid3<-c(1,1,1,3,3,3)
pid3<-c(1,2,3,1,2,3)
bin<-c(1,2,3,1,2,3)
x.breaks<-list(breaks.dems, breaks.dems, breaks.dems, breaks.reps, breaks.reps, breaks.reps)
x.ranges<-list(range.dems, range.dems, range.dems, range.reps, range.reps, range.reps)
titles<-c("Democratic Voter/Democratic Leg.","Independent Voter/Democratic Leg.","Republican Voter/Democratic Leg.","Democratic Voter/Republican Leg.","Independent Voter/Republican Leg.","Republican Voter/Republican Leg.")

graphs<-list()
for (j in 1:6){
	graphs[[j]] <- ggplot(merged2[merged2$pid3==pid3[j] & merged2$leg_pid3==leg_pid3[j] & !is.na(merged2$weight),],  aes( x = iv, y = dv)) + theme (axis.title.y  = element_text(size= rel(1)))
    graphs[[j]] <-graphs[[j]] +ylab(x.text)+ xlab("Incumbent Position") + theme(legend.position="none")+ scale_x_continuous(limits=x.ranges[[j]], breaks=x.breaks[[j]])+scale_y_continuous( breaks=c(0,.2,.4,.6,.8,1))+ ggtitle(titles[j])  +coord_cartesian(ylim=c(0,1)) + stat_smooth(se=TRUE, method="gam", formula = y ~ s(x, k=4), aes(weight=weight), data=merged2[merged2$leg_pid3==leg_pid3[j] & merged2$pid3==pid3[j] & !is.na(merged2$weight),],colour="black",linetype="solid", size=.8)+ theme_bw()+ theme(plot.title = element_text(hjust = 0.5))


	for (p in bin){
		graphs[[j]]  <- graphs[[j]]  + stat_smooth(se=TRUE, method="gam", formula = y ~ s(x, k=4), aes(weight=weight), data=merged2[merged2$bin==bin[p] & merged2$leg_pid3==leg_pid3[j] & merged2$pid3==pid3[j] & !is.na(merged2$weight),],colour="black",linetype=line.type[p], size=.8)+ theme(plot.title = element_text(size = rel(1),hjust = 0.5))
		}
	print("pid group:")
	print(j)
}

pdf(paste(analysis,"_loess_weighted_",graph_types[i],"_simple.pdf",sep=""),width=7,height=4)
multiplot(allvotes_demleg,allvotes_repleg, cols=2)
dev.off()         
	
pdf(paste("Figure1",analysis,"_loess_weighted_",graph_types[i], ".pdf",sep=""),width=7,height=8)
multiplot(graphs[[1]], graphs[[2]], graphs[[3]], graphs[[4]],graphs[[5]],graphs[[6]], cols=2)
dev.off()  



#######################################
### Make regression tables
### (Main Regression in Spatial Voting Paper (January 2015)
#######################################


model_types<-c("all_races" )
titles<-c("Voting for Incumbents")
x.text<-"Pr(Vote for Incumbent)"

merged$congress_post_vote[merged$congress_post_vote>1]<-NA

merged2<-merged
merged2<-merged2[merged2$incumb %in% c(-1,1) & !is.na(merged2$challengers.total.receipts),]


merged2<-subset(merged2, !is.na(pid3) &!is.na(leg_pid3) & !is.na(congress_post_vote) & !is.na(survey_ideology) & !is.na(leg_pid3) & !is.na(leg_ideology) & !is.na(year) & !is.na(cd) & leg_pid3!=2 & year>2004 & !is.na(weight))
colnames(merged2)[colnames(merged2)=="fipsyear_vote.x"]<-"fipsyear_vote"


## CW, September 2016
merged2$rhousevote<-NA
merged2$rhousevote[merged2$congress_post_vote==0]<-1
merged2$rhousevote[merged2$congress_post_vote==1]<-0

merged2$std_survey_ideology<-scale(merged2$survey_ideology)[,1]
merged2$std_leg_ideology<-scale(merged2$leg_ideology)[,1]
merged2$std_survey_ideology_sq<-merged2$std_survey_ideology^2
merged2$std_leg_ideology_sq<-merged2$std_leg_ideology^2

merged2$std_leg_ideology_dems<-NA
merged2$std_leg_ideology_reps<-NA
merged2$std_leg_ideology_sq_dems<-NA
merged2$std_leg_ideology_sq_reps<-NA

merged2$std_leg_ideology_dems[merged2$leg_pid3==1]<-merged2$std_leg_ideology[merged2$leg_pid3==1]
merged2$std_leg_ideology_reps[merged2$leg_pid3==3]<-merged2$std_leg_ideology[merged2$leg_pid3==3]

merged2$std_leg_ideology_sq_dems[merged2$leg_pid3==1]<-merged2$std_leg_ideology_sq[merged2$leg_pid3==1]
merged2$std_leg_ideology_sq_reps[merged2$leg_pid3==3]<-merged2$std_leg_ideology_sq[merged2$leg_pid3==3]

save("merged2",file="data_for_voteshares.RData")

#####
### Table 2, column 1
#####
mod1d<-(lm(rhousevote~std_survey_ideology+std_survey_ideology_sq+std_leg_ideology_dems+std_leg_ideology_sq_dems+std_survey_ideology*std_leg_ideology_dems+factor(pid3)+factor(source), data=merged2[merged2 $leg_pid3==1,], weight=weight))
summary(mod1d)
mod1d.cluster=coeftest(mod1d,vcov=vcovCluster(mod1d,merged2$fipsyear_vote[!is.na(merged2$std_survey_ideology)  & !is.na(merged2$std_leg_ideology_dems) & !is.na(merged2$rhousevote) & !is.na(merged2$weight)]))

#####
### Table 2, column 2
#####
mod1r<-(lm(rhousevote~std_survey_ideology+std_survey_ideology_sq+std_leg_ideology_reps+std_leg_ideology_sq_reps+std_survey_ideology*std_leg_ideology_reps+factor(pid3)+factor(source), data=merged2[merged2 $leg_pid3==3,], weight=weight))
summary(mod1r)
mod1r.cluster=coeftest(mod1r,vcov=vcovCluster(mod1r,merged2$fipsyear_vote[!is.na(merged2$std_survey_ideology)  & !is.na(merged2$std_leg_ideology_reps) & !is.na(merged2$rhousevote) & !is.na(merged2$weight)]))

### replicate analysis with DW Dime Measures

## dwdime already standard, normal so don't need to standardize
merged2$std_leg_ideology_dwdime_dems<-(merged2$dwdime_dem)
merged2$std_leg_ideology_dwdime_reps<-(merged2$dwdime_rep)
merged2$std_leg_ideology_dwdime_sq_dems<-merged2$std_leg_ideology_dwdime_dems^2
merged2$std_leg_ideology_dwdime_sq_reps<-merged2$std_leg_ideology_dwdime_reps^2


#####
### Table 2, column 3
#####
mod3<- lm(rhousevote ~ std_survey_ideology+std_survey_ideology_sq+ std_leg_ideology_dwdime_dems+std_leg_ideology_dwdime_sq_dems +std_leg_ideology_dwdime_reps+std_leg_ideology_dwdime_sq_reps +std_survey_ideology*std_leg_ideology_dwdime_dems+
std_survey_ideology*std_leg_ideology_dwdime_reps + factor(pid3)+ factor(source)  ,data=merged2[], weight=weight)
summary(mod3)
mod3.cluster=coeftest(mod3,vcov=vcovCluster(mod3,merged2$fipsyear_vote[!is.na(merged2$std_survey_ideology)  & !is.na(merged2$std_leg_ideology_dwdime_sq_dems)& !is.na(merged2$std_leg_ideology_dwdime_sq_reps) & !is.na(merged2$rhousevote) & !is.na(merged2$weight)]))



#####
### Table 2, column 4
#####
### replicate Adams et al with our ideology measure
data$source<-"CCES_2010"


data2<-merge(data, merged2, by.x=c("cces_id", "source"),by.y=c("respondent", "source"))

data2$std_leg_ideology_dems<-scale(data2$d_cand_ip)[,1]
data2$std_leg_ideology_reps<-scale(data2$r_cand_ip)[,1]
data2$std_leg_ideology_sq_dems<-data2$std_leg_ideology_dems^2
data2$std_leg_ideology_sq_reps<-data2$std_leg_ideology_reps^2

normal_range <- data2$d_cand_ip < 0.2 & data2$r_cand_ip > 0.3 &
                data2$r_cand_ip < 1.62

data2$pid3<-data2$pid3.y

mod2 <- lm(rhousevote.x ~ std_survey_ideology+std_survey_ideology_sq+ std_leg_ideology_dems+std_leg_ideology_sq_dems +std_leg_ideology_reps+std_leg_ideology_sq_reps +std_survey_ideology*std_leg_ideology_dems+
std_survey_ideology*std_leg_ideology_reps + factor(pid3) ,data=data2[], weight=weight.x)
summary(mod2)
mod2.cluster=coeftest(mod2,vcov=vcovCluster(mod2,data2$fipsyear_vote[!is.na(data2$std_survey_ideology) &!is.na(data2$std_leg_ideology_dems) & !is.na(data2$std_leg_ideology_reps) & !is.na(data2$rhousevote.x) & !is.na(data2$weight.x)]))


### MAIN RESULTS in TABLE 4
stargazer(mod1d.cluster,mod1r.cluster, mod3.cluster, mod2.cluster)
stargazer(mod1d,mod1r, mod3, mod2)



##############################################	
## Main graph with substantive effects of legislator ideology 
## (Figures 2 and 3)
##############################################	

## Change working directory to RData folder so everything saves to the right place


### Predicted Values
colnames(merged2)[colnames(merged2)=="fipsyear_vote.x"]<-"fipsyear_vote"

merged2$leg_dem<-NA
merged2$leg_dem[merged2$leg_pid3==1]<-1
merged2$leg_dem[merged2$leg_pid3==3]<-0

merged2$congress_post_vote[merged2$congress_post_vote>1]<-NA
merged2<-merged2[merged2$incumb %in% c(-1,1) & !is.na(merged2$challengers.total.receipts),]
merged2$dv<-merged2$ivote_house


#########
## The following pieces of code estimate regression and boostrapped models for generating standard errors of predicted values 
## errors for a series of models. Then it graphs tehm.
#########


#########
## Predicted values for main model
#########

		
sims<- 20
model.results.dem.legs<-vector(mode="list", length=sims)
model.results.rep.legs<-vector(mode="list", length=sims)
model.results.bothcands<-vector(mode="list", length=sims)
model.results.linear.left<-vector(mode="list", length=sims)
model.results.linear.middle<-vector(mode="list", length=sims)
model.results.linear.right<-vector(mode="list", length=sims)
model.results.bothcands.dwdime<-vector(mode="list", length=sims)


model.predictions.dems<-matrix(nrow=0, ncol=3)
model.predictions.reps<-matrix(nrow=0, ncol=3)
colnames(model.predictions.dems)<-c("ideo", "prediction", "pid3")
colnames(model.predictions.reps)<-c("ideo", "prediction", "pid3")

model.predictions.dems.both<-matrix(nrow=0, ncol=3)
model.predictions.reps.both<-matrix(nrow=0, ncol=3)
colnames(model.predictions.dems.both)<-c("ideo", "prediction", "pid3")
colnames(model.predictions.reps.both)<-c("ideo", "prediction", "pid3")

model.predictions.dems.both.dwdime<-matrix(nrow=0, ncol=3)
model.predictions.reps.both.dwdime<-matrix(nrow=0, ncol=3)
colnames(model.predictions.dems.both.dwdime)<-c("ideo", "prediction", "pid3")
colnames(model.predictions.reps.both.dwdime)<-c("ideo", "prediction", "pid3")

model.predictions.dems.linear<-matrix(nrow=0, ncol=3)
model.predictions.reps.linear<-matrix(nrow=0, ncol=3)
colnames(model.predictions.dems.linear)<-c("ideo", "prediction", "pid3")
colnames(model.predictions.reps.linear)<-c("ideo", "prediction", "pid3")

data2$d_cand_ip_sq<-data2$d_cand_ip^2
data2$r_cand_ip_sq<-data2$r_cand_ip^2
for (i in 1:sims){
		boot2<-data2[sample(1:dim(data2)[1], dim(data2)[1], replace=T),]
		model.bothcandidates<-glm(rhousevote.x ~ std_survey_ideology* factor(pid3)+std_survey_ideology_sq* factor(pid3)+ d_cand_ip* factor(pid3)+d_cand_ip_sq* factor(pid3) +r_cand_ip* factor(pid3)+r_cand_ip_sq* factor(pid3) +std_survey_ideology*d_cand_ip* factor(pid3)+
std_survey_ideology*r_cand_ip* factor(pid3),data=boot2[],weight=weight.x,family=binomial(link = "logit"))
		model.results.bothcands[[i]]<-model.bothcandidates

		boot2<-merged2[sample(1:dim(merged2)[1], dim(merged2)[1], replace=T),]
		model.bothcandidates.dwdime<-glm(rhousevote ~ std_survey_ideology+std_survey_ideology_sq+ std_leg_ideology_dwdime_dems+std_leg_ideology_dwdime_sq_dems +std_leg_ideology_dwdime_reps+std_leg_ideology_dwdime_sq_reps +std_survey_ideology*std_leg_ideology_dwdime_dems+
std_survey_ideology*std_leg_ideology_dwdime_reps + factor(pid3) ,data=boot2[], weight=weight,family=binomial(link = "logit"))
		model.results.bothcands.dwdime[[i]]<-model.bothcandidates.dwdime
	
	for (p in 0:1){
		print(paste("i: ", i, sep=""))
		boot<-merged2[sample(1:dim(merged2)[1], dim(merged2)[1], replace=T),]
		model.incumbents<-(glm(rhousevote~std_survey_ideology*factor(pid3)+ std_survey_ideology_sq*factor(pid3)+std_leg_ideology*factor(pid3)+std_leg_ideology_sq*factor(pid3)+std_survey_ideology*std_leg_ideology*factor(pid3), data=boot[boot $leg_dem==p,],weight=weight,family=binomial(link = "logit")))	

		if (p==0){
		model.results.rep.legs[[i]]<-model.incumbents
		}
		if (p==1){
		model.results.dem.legs[[i]]<-model.incumbents
		}
	}
}

### All survey data, quadratic model
survey_ideo<-c(-1,0,1)
for (p in 0:1){
	if(p==0){
		for (i in 1:sims){
			for (j in 1:3){
				temp<-NA
				temp$ideo<-seq(-1.7,-.1, .1)
				temp$prediction<-predict(model.results.dem.legs[[i]], data.frame(std_leg_ideology =seq(-1.7,-.1, .1), std_leg_ideology_sq = seq(-1.7,-.1, .1)^2, std_survey_ideology=survey_ideo[j],std_survey_ideology_sq=survey_ideo[j]^2,  pid3=j), type = "response")
				temp$pid3<-j
				temp<-as.data.frame(temp)
				model.predictions.dems<-rbind(model.predictions.dems, temp[,2:4])
			}
		}
	}

	if(p==1){
		for (i in 1:sims){
			for (j in 1:3){
				temp<-NA
				temp$ideo<-seq(.1,1.7, .1)
				temp$prediction<-predict(model.results.rep.legs[[i]], data.frame(std_leg_ideology =seq(.1,1.7, .1), std_leg_ideology_sq = seq(.1,1.7, .1)^2, std_survey_ideology=survey_ideo[j],std_survey_ideology_sq=survey_ideo[j]^2,  pid3=j), type = "response")
				temp$pid3<-j
				temp<-as.data.frame(temp)
				model.predictions.reps<-rbind(model.predictions.reps, temp[,2:4])
			}
		}
	}
}


## Adams et al data, quadratic model
survey_ideo<-c(-1,0,1)
for (p in 0:1){
	if(p==0){
		for (i in 1:sims){
			for (j in 1:3){
				temp<-NA
				temp$ideo<-seq(-1.7,-.1, .1)
				temp$prediction<-predict(model.results.bothcands[[i]], data.frame(d_cand_ip =seq(-1.7,-.1, .1), d_cand_ip_sq = seq(-1.7,-.1, .1)^2, r_cand_ip=1, r_cand_ip_sq=1, std_survey_ideology=survey_ideo[j],std_survey_ideology_sq=survey_ideo[j]^2,  pid3=j), type = "response")
				temp$pid3<-j
				temp<-as.data.frame(temp)
				model.predictions.dems.both<-rbind(model.predictions.dems.both, temp[,2:4])
			}
		}
	}

	if(p==1){
		for (i in 1:sims){
			for (j in 1:3){
				temp<-NA
				temp$ideo<-seq(.1,1.7, .1)
				temp$prediction<-predict(model.results.bothcands[[i]], data.frame(r_cand_ip =seq(.1,1.7, .1), r_cand_ip_sq = seq(.1,1.7, .1)^2, d_cand_ip=-1, d_cand_ip_sq=1, std_survey_ideology=survey_ideo[j],std_survey_ideology_sq=survey_ideo[j]^2,  pid3=j), type = "response")
				temp$pid3<-j
				temp<-as.data.frame(temp)
				model.predictions.reps.both<-rbind(model.predictions.reps.both, temp[,2:4])
				print(j)
			}
		}
	}
}


## Bonica DW-Dime data, quadratic model
survey_ideo<-c(-1,0,1)
for (p in 0:1){
	if(p==0){
		for (i in 1:sims){
			for (j in 1:3){
				temp<-NA
				temp$ideo<-seq(-1.5,-.1, .1)
				temp$prediction<-predict(model.results.bothcands.dwdime[[i]], data.frame( std_leg_ideology_dwdime_dems =seq(-1.5,-.1, .1), std_leg_ideology_dwdime_sq_dems = seq(-1.5,-.1, .1)^2,std_leg_ideology_dwdime_reps =1, std_leg_ideology_dwdime_sq_reps =1^2  ,std_survey_ideology=survey_ideo[j],std_survey_ideology_sq=survey_ideo[j]^2,  pid3=j), type = "response")
				temp$pid3<-j
				temp<-as.data.frame(temp)
				model.predictions.dems.both.dwdime<-rbind(model.predictions.dems.both.dwdime, temp[,2:4])
			}
		}
	}

	if(p==1){
		for (i in 1:sims){
			for (j in 1:3){
				temp<-NA
				temp$ideo<-seq(.1,1.7, .1)
				temp$prediction<-predict(model.results.bothcands.dwdime[[i]], data.frame( std_leg_ideology_dwdime_dems =-.86, std_leg_ideology_dwdime_sq_dems =(-.86)^2,std_leg_ideology_dwdime_reps =seq(.1,1.7, .1), std_leg_ideology_dwdime_sq_reps =seq(.1,1.7, .1)^2,  std_survey_ideology=survey_ideo[j],std_survey_ideology_sq=survey_ideo[j]^2,  pid3=j), type = "response")
				temp$pid3<-j
				temp<-as.data.frame(temp)
				model.predictions.reps.both.dwdime<-rbind(model.predictions.reps.both.dwdime, temp[,2:4])
				print(j)
			}
		}
	}
}



pdf(file="Figure2b_spatial_voting_rep_incumbent.pdf")
reps.graph <- ggplot(model.predictions.reps, aes(ideo, prediction,  colour=factor(pid3),  
					linetype = factor(pid3)))
reps.graph <- reps.graph+ stat_smooth()
reps.graph <- reps.graph + ylab('Pr(Vote for Republican)')
reps.graph <- reps.graph + xlab('Republican Legislator Ideology')
#reps.graph <- reps.graph +  scale_colour_manual(values = c("blue","green", "red"))
reps.graph<- reps.graph+ theme_bw()
#reps.graph <- reps.graph + theme( legend.position="bottom")
reps.graph <- reps.graph +   scale_linetype_manual(values=c("dotdash", "longdash", "dashed"))
reps.graph <- reps.graph + theme(legend.position="none")
reps.graph <- reps.graph +scale_colour_manual(values=c("black","black", "black"))
#						name="Voter Party ID",
#                       breaks=c("1", "2", "3"),
 #                      labels=c("Democrat", "Independent", "Republican"))
   reps.graph <- reps.graph +  theme(axis.title.x = element_text(size=16))
reps.graph <- reps.graph +  theme(axis.title.y = element_text(size=16))
 #  reps.graph <- reps.graph + theme(legend.text=element_text(size=13))
reps.graph <- reps.graph +  scale_y_continuous(limits = c(-0.05, 1.05), breaks=seq(0,1,.1))
reps.graph <- reps.graph +   theme(text = element_text(size=16))
reps.graph
dev.off()


pdf(file="Figure2a_spatial_voting_dem_incumbent.pdf")
dems.graph <- ggplot(model.predictions.dems, aes(ideo, prediction,  colour=factor(pid3),  
					linetype = factor(pid3)))
dems.graph <- dems.graph+ stat_smooth()
dems.graph <- dems.graph + ylab('Pr(Vote for Republican)')
dems.graph <- dems.graph + xlab('Democratic Legislator Ideology')
#dems.graph <- dems.graph +  scale_colour_manual(values = c("blue","green", "red"))
dems.graph<- dems.graph+ theme_bw()
dems.graph <- dems.graph +   scale_linetype_manual(values=c("dotdash", "longdash", "dashed"))
#dems.graph <- dems.graph + theme( legend.position="bottom")
dems.graph <- dems.graph + theme(legend.position="none")
dems.graph <- dems.graph +scale_colour_manual(values=c("black","black", "black"))
#						name="Voter Party ID",
 #                      breaks=c("1", "2", "3"),
 #                      labels=c("Democrat", "Independent", "Republican"))
#dems.graph <-  dems.graph +scale_linetype_manual(values = c('solid', 'longdash', 'dashed'),
#						name="Voter Party ID",
 #                      labels=c("Democrat", "Independent", "Republican")) 
   dems.graph <- dems.graph +  theme(axis.title.x = element_text(size=16))
   dems.graph <- dems.graph +  theme(axis.title.y = element_text(size=16))
dems.graph <- dems.graph +  scale_y_continuous(limits = c(-0.05, 1.05), breaks=seq(0,1,.1))
   #dems.graph <- dems.graph + theme(legend.text=element_text(size=13))
dems.graph <- dems.graph +   theme(text = element_text(size=16))
dems.graph <-dems.graph 
dems.graph
dev.off()

## Adams et al data
pdf(file="spatial_voting_dem_candidates_adamsetal.pdf")
dems.graph <- ggplot(model.predictions.dems.both, aes(ideo, prediction,  colour=factor(pid3)))
dems.graph <- dems.graph+ stat_smooth()
dems.graph <- dems.graph + ylab('Pr(Vote for Republican)')
dems.graph <- dems.graph + xlab('Democratic Candidate Ideology')
#dems.graph <- dems.graph +  scale_colour_manual(values = c("blue","green", "red"))
dems.graph<- dems.graph+ theme_bw()
dems.graph <- dems.graph + theme( legend.position="bottom")
dems.graph <- dems.graph +scale_colour_manual(values=c("blue","green", "red"),  
						name="Voter Party ID",
                       breaks=c("1", "2", "3"),
                       labels=c("Democrat", "Independent", "Republican"))
dems.graph <- dems.graph +  scale_y_continuous(limits = c(-0.05, 1.05), breaks=seq(0,1,.1))
   dems.graph <- dems.graph +  theme(axis.title.x = element_text(size=14))
   dems.graph <- dems.graph +  theme(axis.title.y = element_text(size=14))
   dems.graph <- dems.graph + theme(legend.text=element_text(size=13))
dems.graph <-dems.graph 
dems.graph
dev.off()

pdf(file="spatial_voting_rep_candidates_adamsetal.pdf")
reps.graph <- ggplot(model.predictions.reps.both, aes(ideo, prediction,  colour=factor(pid3)))
reps.graph <- reps.graph+ stat_smooth()
reps.graph <- reps.graph + ylab('Pr(Vote for Republican)')
reps.graph <- reps.graph + xlab('Republican Candidate Ideology')
#reps.graph <- reps.graph +  scale_colour_manual(values = c("blue","green", "red"))
reps.graph<- reps.graph+ theme_bw()
reps.graph <- reps.graph + theme( legend.position="bottom")
reps.graph <- reps.graph +scale_colour_manual(values=c("blue","green", "red"),  
						name="Voter Party ID",
                       breaks=c("1", "2", "3"),
                       labels=c("Democrat", "Independent", "Republican"))
reps.graph <- reps.graph +  scale_y_continuous(limits = c(-0.05, 1.05), breaks=seq(0,1,.1))
   reps.graph <- reps.graph +  theme(axis.title.x = element_text(size=14))
   reps.graph <- reps.graph +  theme(axis.title.y = element_text(size=14))
   reps.graph <- reps.graph + theme(legend.text=element_text(size=13))
reps.graph <-reps.graph +   theme(plot.title = element_text(size = 10),axis.title.y =
                     element_text(size = 14), axis.text =
                 element_text(size = 10), axis.title.x =
                     element_text(size = 14))
reps.graph
dev.off()

## DW-Dime data

pdf(file="Figure3a_spatial_voting_dem_candidates_dwdime.pdf")
dems.graph <- ggplot(model.predictions.dems.both.dwdime, aes(ideo, prediction,  colour=factor(pid3),  
					linetype = factor(pid3)))
dems.graph <- dems.graph+ stat_smooth()
dems.graph <- dems.graph + ylab('Pr(Vote for Republican)')
dems.graph <- dems.graph + xlab('Democratic Candidate Ideology')
#dems.graph <- dems.graph +  scale_colour_manual(values = c("blue","green", "red"))
dems.graph<- dems.graph+ theme_bw()
#dems.graph <- dems.graph + theme( legend.position="bottom")
dems.graph <- dems.graph + theme(legend.position="none")
dems.graph <- dems.graph +   scale_linetype_manual(values=c("dotdash", "longdash", "dashed"))
dems.graph <- dems.graph +scale_colour_manual(values=c("black","black", "black"))  
	#					name="Voter Party ID",
    #                   breaks=c("1", "2", "3"),
   #                    labels=c("Democrat", "Independent", "Republican"))
dems.graph <- dems.graph +  scale_y_continuous(limits = c(-0.05, 1.05), breaks=seq(0,1,.1))
   dems.graph <- dems.graph +  theme(axis.title.x = element_text(size=16))
   dems.graph <- dems.graph +  theme(axis.title.y = element_text(size=16))
 #  dems.graph <- dems.graph + theme(legend.text=element_text(size=13))
dems.graph <-dems.graph 
dems.graph <- dems.graph +   theme(text = element_text(size=16))
dems.graph
dev.off()

pdf(file="Figure3b_spatial_voting_rep_candidates_dwdime.pdf")
reps.graph <- ggplot(model.predictions.reps.both.dwdime, aes(ideo, prediction,  colour=factor(pid3),  
					linetype = factor(pid3)))
reps.graph <- reps.graph+ stat_smooth()
reps.graph <- reps.graph + ylab('Pr(Vote for Republican)')
reps.graph <- reps.graph + xlab('Republican Candidate Ideology')
reps.graph <- reps.graph +   scale_linetype_manual(values=c("dotdash", "longdash", "dashed"))
#reps.graph <- reps.graph +  scale_colour_manual(values = c("blue","green", "red"))
reps.graph<- reps.graph+ theme_bw()
#reps.graph <- reps.graph + theme( legend.position="bottom")
reps.graph <- reps.graph + theme(legend.position="none")

reps.graph <- reps.graph +scale_colour_manual(values=c("black","black", "black"))  
#						name="Voter Party ID",
 #                      breaks=c("1", "2", "3"),
  #                     labels=c("Democrat", "Independent", "Republican"))
reps.graph <- reps.graph +  scale_y_continuous(limits = c(-0.05, 1.05), breaks=seq(0,1,.1))
   reps.graph <- reps.graph +  theme(axis.title.x = element_text(size=16))
   reps.graph <- reps.graph +  theme(axis.title.y = element_text(size=16))
  # reps.graph <- reps.graph + theme(legend.text=element_text(size=13))
reps.graph <- reps.graph +   theme(text = element_text(size=16))
reps.graph
dev.off()




### MAIN RESULTS in TABLE 5


## joint scaling model

data$left <- ifelse(data$np_score < data$d_cand_ip ,1,0)
                      
data$right <- ifelse(data$np_score > data$r_cand_ip,1,0)
                                            
data$middle <- ifelse(data$np_score > data$d_cand_ip &
                      data$np_score < data$r_cand_ip,1,0)

mod.left <- lm(rhousevote ~ np_score +r_cand_ip + d_cand_ip + factor(pid3) + repseat + incumbent_party + attendchurch + female + black + latino + age + college_grad + own_home,data=data[data$left==1,], weight=weight) # rep_spending_advantage

mod.left.cluster=coeftest(mod.left,vcov=vcovCluster(mod.left,data$district[!is.na(data$np_score) &!is.na(data$r_cand_ip) & !is.na(data$d_cand_ip)& !is.na(data$pid3)  & !is.na(data$rhousevote) & !is.na(data$weight)& !is.na(data$repseat)& !is.na(data$incumbent_party)& !is.na(data$attendchurch)& !is.na(data$female)& !is.na(data$black) & !is.na(data$latino) & !is.na(data$age) & !is.na(data$college_grad) & !is.na(data$own_home)& data$left==1]))

mod.middle <- lm(rhousevote ~ np_score +r_cand_ip + d_cand_ip + factor(pid3) + repseat + incumbent_party + attendchurch + female + black + latino + age + college_grad + own_home,data=data[data$middle==1,], weight=weight) # rep_spending_advantage

mod.middle.cluster=coeftest(mod.middle,vcov=vcovCluster(mod.middle,data$district[!is.na(data$np_score) &!is.na(data$r_cand_ip) & !is.na(data$d_cand_ip)& !is.na(data$pid3)  & !is.na(data$rhousevote) & !is.na(data$weight)& !is.na(data$repseat)& !is.na(data$incumbent_party)& !is.na(data$attendchurch)& !is.na(data$female)& !is.na(data$black) & !is.na(data$latino) & !is.na(data$age) & !is.na(data$college_grad) & !is.na(data$own_home)& data$middle==1]))

mod.right <- lm(rhousevote ~ np_score +r_cand_ip + d_cand_ip + factor(pid3) + repseat + incumbent_party + attendchurch + female + black + latino + age + college_grad + own_home,data=data[data$right==1,], weight=weight) # rep_spending_advantage

mod.right.cluster=coeftest(mod.right,vcov=vcovCluster(mod.right,data$district[!is.na(data$np_score) &!is.na(data$r_cand_ip) & !is.na(data$d_cand_ip)& !is.na(data$pid3)  & !is.na(data$rhousevote) & !is.na(data$weight)& !is.na(data$repseat)& !is.na(data$incumbent_party)& !is.na(data$attendchurch)& !is.na(data$female)& !is.na(data$black) & !is.na(data$latino) & !is.na(data$age) & !is.na(data$college_grad) & !is.na(data$own_home)& data$right==1]))

stargazer(mod.left,mod.middle, mod.right)
stargazer(mod.left.cluster,mod.middle.cluster, mod.right.cluster)

mod.left <- glm(rhousevote ~ np_score +r_cand_ip + d_cand_ip + factor(pid3) + repseat + incumbent_party + attendchurch + female + black + latino + age + college_grad + own_home,family=binomial(link="logit"),data=data[data$left==1,], weight=weight) # rep_spending_advantage
mod.middle <- glm(rhousevote ~ np_score +r_cand_ip + d_cand_ip + factor(pid3) + repseat + incumbent_party + attendchurch + female + black + latino + age + college_grad + own_home,family=binomial(link="logit"),data=data[data$middle==1,], weight=weight) # rep_spending_advantage
mod.right <- glm(rhousevote ~ np_score +r_cand_ip + d_cand_ip + factor(pid3) + repseat + incumbent_party + attendchurch + female + black + latino + age + college_grad + own_home,family=binomial(link="logit"),data=data[data$right==1,], weight=weight) # rep_spending_advantage
#stargazer(mod.left,mod.middle, mod.right)
   
mod <- glm(rhousevote ~ np_score:left + np_score:middle +np_score:right+ r_cand_ip*right + d_cand_ip*left + pid7 + repseat + incumbent_party + attendchurch + female + black + latino + age + college_grad + own_home,family=binomial(link="logit"),data=data, weight=weight) # rep_spending_advantage
summary(mod)

#####
### Republican candidate moves, and Democratic candidate stays fixed
#####

## voter to the left of Democrat
mod <- glm(rhousevote ~ np_score + r_cand_ip + d_cand_ip + pid3,family=binomial(link="logit"),data=data[data$left==1,], weight=weight) # rep_spending_advantage
prediction<-predict.glm(mod, data.frame(np_score = -1.075,
   r_cand_ip = seq(.4,1.6,.1), d_cand_ip = -0.660,left=1, middle=0, right=0, pid3=-1),
   type = "response")
prediction   
## voter between Democrat and Republican
mod <- glm(rhousevote ~ np_score + r_cand_ip + d_cand_ip + pid3,family=binomial(link="logit"),data=data[data$middle==1,], weight=weight) # rep_spending_advantage
prediction<-predict.glm(mod, data.frame(np_score = 0.134 ,
   r_cand_ip = seq(.4,1.6,.1), d_cand_ip = -0.660,left=0, middle=1, right=0, pid3=0),
   type = "response")
prediction   
## voter to the right of Republican
mod <- glm(rhousevote ~ np_score + r_cand_ip + d_cand_ip + pid3,family=binomial(link="logit"),data=data[data$right==1,], weight=weight) # rep_spending_advantage
prediction<-predict.glm(mod, data.frame(np_score = 1.228 ,
   r_cand_ip = seq(.4,1.6,.1), d_cand_ip = -0.660,left=0, middle=0, right=1, pid3=1),
   type = "response")
prediction   

#####
### Democratic candidate moves, and Republican candidate stays fixed
#####

## voter to the left of Democrat
mod <- glm(rhousevote ~ np_score + r_cand_ip + d_cand_ip + pid3,family=binomial(link="logit"),data=data[data$left==1,], weight=weight) # rep_spending_advantage
prediction<-predict.glm(mod, data.frame(np_score = -1.075,
   d_cand_ip = seq(-1.3,.13,.1), r_cand_ip = 1.066,left=1, middle=0, right=0, pid3=-1),
   type = "response")
prediction   
## voter between Democrat and Republican
mod <- glm(rhousevote ~ np_score + r_cand_ip + d_cand_ip + pid3,family=binomial(link="logit"),data=data[data$middle==1,], weight=weight) # rep_spending_advantage
prediction<-predict.glm(mod, data.frame(np_score = 0.134 ,
   d_cand_ip = seq(-1.3,.13,.1), r_cand_ip = 1.066,left=0, middle=1, right=0, pid3=0),
   type = "response")
prediction   
## voter to the right of Republican
mod <- glm(rhousevote ~ np_score + r_cand_ip + d_cand_ip + pid3,family=binomial(link="logit"),data=data[data$right==1,], weight=weight) # rep_spending_advantage
prediction<-predict.glm(mod, data.frame(np_score = 1.228 ,
   d_cand_ip = seq(-1.3,.13,.1), r_cand_ip = 1.066,left=0, middle=0, right=1, pid3=1),
   type = "response")
prediction   

stargazer(mod)



##############################################	
## Generate Predicted Vote Shares 
## (Figure 4)
##############################################	

###############
load("dem_vote_prediction.RData")

above_to_below<-0

mat<-subset(mat, !is.na(mat[,1]))
mat<-as.data.frame(mat)

samp<-sample(dim(mat)[1], 50, replace = FALSE, prob = NULL)

mat2<-mat[samp,]
mat3<-matrix(data=NA, nrow=0, ncol=3)
colnames(mat3)<-c("leg_id", "vote_share", "position")
mat3<-as.data.frame(mat3)
for (i in 1:dim(mat2)[1]){
	temp<-as.numeric(as.vector(mat2[i,]))
	temp2<-cbind(rownames(mat2)[i],temp, colnames(mat2))
	mat3<-rbind(mat3, temp2)
}
colnames(mat3)<-c("leg_id", "vote_share", "position")
leg<-unique(mat3$leg_id)
mat3$position<-as.numeric(as.vector(mat3$position))
mat3$vote_share<-as.numeric(as.vector(mat3$vote_share))

one_sd_change<-as.numeric(mat[,8])-as.numeric(mat[,5])
one_sd_change<-as.data.frame(one_sd_change)
mean(as.numeric(as.vector(one_sd_change[,1])))

x.axis<-as.numeric(colnames(mat))
grab<-1:dim(mat)[1]
a <- ggplot(mat3)
a <-a + ylab("Projected Vote Share")+ xlab("Incumbent Position") + theme(legend.position="none")+theme(axis.title=element_text(size=10))+scale_y_continuous(limits = c(.3,.9))+ ggtitle("Democratic Incumbents")  +theme(title=element_text(size=9))

	for (p in 1:length(leg)){
		a <- a + geom_smooth(data=mat3[mat3$leg_id==leg[p],], aes(x=position, y = vote_share), colour="black", size = .25, se = FALSE, stat = "smooth")
		print(p)
	}
a <-a + geom_hline(aes(yintercept=.5),size = 1.5, colour="black")+ theme_bw()+theme(plot.title=element_text(size=10))+theme(axis.title=element_text(size=10))

theme_bare <- theme(
  axis.line = element_blank(), 
  axis.text.y = element_blank(),
  axis.title.y = element_blank(), 
  panel.grid.major = element_blank(),
  panel.grid.minor = element_blank()
)


b <- ggplot(one_sd_change, aes(x = one_sd_change))
b <- b + geom_density(colour="black", fill="grey")+ xlab("Democrats' Vote Share Gain") + ylab("")+theme_bare+ theme_bw()+theme(plot.title=element_text(size=10))+theme(axis.title=element_text(size=10))

one_sd_change2<-as.data.frame(leg_vs_sd_mod-leg_vs)
colnames(one_sd_change2)<-"one_sd_change"
mean(as.numeric(as.vector(one_sd_change2[,1])))

b2 <- ggplot(one_sd_change2, aes(x = one_sd_change))
b2 <- b2 + geom_density(colour="black", fill="grey")+ xlab("Democrats' Vote Share Gain") + ylab("")+theme_bare+ theme_bw()+theme(plot.title=element_text(size=10))+theme(axis.title=element_text(size=10))+ theme( axis.text.y = element_blank())

d2 <- ggplot(one_sd_change2, aes(x = one_sd_change))
d2 <- d2 + geom_density(colour="black", fill="grey")+ xlab("Republicans' Vote Share Gain") + ylab("")+theme_bare+ theme_bw()+theme(plot.title=element_text(size=10))+theme(axis.title=element_text(size=10))+ theme(
  axis.text.y = element_blank())


load("rep_vote_prediction.RData")

mat<-subset(mat, !is.na(mat[,1]))
mat<-as.data.frame(mat)

one_sd_change<-as.numeric(mat[,5])-as.numeric(mat[,8])
one_sd_change<-as.data.frame(one_sd_change)
mean(as.numeric(as.vector(one_sd_change[,1])))

samp<-sample(dim(mat)[1], 50, replace = FALSE, prob = NULL)

mat2<-mat[samp,]
mat3<-matrix(data=NA, nrow=0, ncol=3)
colnames(mat3)<-c("leg_id", "vote_share", "position")
mat3<-as.data.frame(mat3)
for (i in 1:dim(mat2)[1]){
	temp<-as.numeric(as.vector(mat2[i,]))
	temp2<-cbind(rownames(mat2)[i],temp, colnames(mat2))
	mat3<-rbind(mat3, temp2)
}
colnames(mat3)<-c("leg_id", "vote_share", "position")
leg<-unique(mat3$leg_id)
mat3$position<-as.numeric(as.vector(mat3$position))
mat3$vote_share<-as.numeric(as.vector(mat3$vote_share))

x.axis<-as.numeric(colnames(mat))
grab<-1:dim(mat)[1]
c<- ggplot(mat3)
c<-c + ylab("Projected Vote Share")+ xlab("Incumbent Position") + theme(legend.position="none")+theme(axis.title=element_text(size=9))+scale_y_continuous(limits = c(.3,.9))+ ggtitle("Republican Incumbents")+ theme_bw()+theme(plot.title=element_text(size=10))+theme(axis.title=element_text(size=10))

	for (p in 1:length(leg)){
		c <- c + geom_smooth(data=mat3[mat3$leg_id==leg[p],], aes(x=position, y = vote_share), colour="black", size = .25, se = FALSE, stat = "smooth")
		print(p)
	}
c <-c + geom_hline(aes(yintercept=.5),size = 1.5, colour="black")

d <- ggplot(one_sd_change, aes(x = one_sd_change))
d <- d + geom_density(colour="black", fill="grey")+ xlab("Republicans' Vote Share Gain") + ylab("")+theme_bare+ theme_bw()+theme(plot.title=element_text(size=10))+theme(axis.title=element_text(size=10))


one_sd_change2<-as.data.frame(leg_vs_sd_mod-leg_vs)
colnames(one_sd_change2)<-"one_sd_change"
mean(as.numeric(as.vector(one_sd_change2[,1])))

one_sd_change2<-as.data.frame(leg_vs_sd_mod-leg_vs)
colnames(one_sd_change2)<-"one_sd_change"
d2 <- ggplot(one_sd_change2, aes(x = one_sd_change))
d2 <- d2 + geom_density(colour="black", fill="grey")+ xlab("Republicans' Vote Share Gain") + ylab("")+theme_bare+ theme_bw()+theme(plot.title=element_text(size=10))+theme(axis.title=element_text(size=10))+ theme(
  axis.text.y = element_blank())

	pdf(paste("Figure4_projected_voteshare_ggplot",".pdf",sep=""),width=8,height=6)
multiplot(b2, d2, cols=2)
dev.off()


##########
### Supplementary Analyses not the paper
##########

## Competitive districts

data2<-merge(data2, cd_pres_results,by.x=c("year","cd_vote"), by.y=c("year","cd_fips"),all.x=T)

mod2.comp <- lm(rhousevote.x ~ std_survey_ideology+std_survey_ideology_sq+ std_leg_ideology_dems+std_leg_ideology_sq_dems +std_leg_ideology_reps+std_leg_ideology_sq_reps +std_survey_ideology*std_leg_ideology_dems+
std_survey_ideology*std_leg_ideology_reps + factor(pid3) ,data=data2[data2$pres_dem >.425 & data2$pres_dem <.575,], weight=weight.x)
summary(mod2.comp)

mod2.comp2 <- lm(rhousevote.x ~ std_survey_ideology+std_survey_ideology_sq+ std_leg_ideology_dems+std_leg_ideology_sq_dems +std_leg_ideology_reps+std_leg_ideology_sq_reps +std_survey_ideology*std_leg_ideology_dems+
std_survey_ideology*std_leg_ideology_reps + factor(pid3) ,data=data2[data2$challengers.total.receipts > 500000,], weight=weight.x)
summary(mod2.comp2)


stargazer(mod2, mod2.comp, mod2.comp2)

## CQ - Key votes 

data2<-merge(data2, cd_pres_results,by.x=c("year","cd_vote"), by.y=c("year","cd_fips"),all.x=T)

mod2.dems.imp <- (lm(rhousevote~std_survey_ideology+ std_survey_ideology_sq+cqkeyvotes_ideal+std_leg_ideology_sq_reps+std_survey_ideology*cqkeyvotes_ideal+factor(pid3)+factor(source), data=merged2[merged2 $leg_pid3==3,], weight=weight))
summary(mod2.dems.imp)

mod2.reps.imp <- lm(rhousevote.x ~ cqkeyvotes_ideal+cqkeyvotes_ideal_sq+ std_leg_ideology_dems+std_leg_ideology_sq_dems +std_leg_ideology_reps+std_leg_ideology_sq_reps +std_survey_ideology*cqkeyvotes_ideal_sq + factor(pid3) ,data=data2[], weight=weight.x)
summary(mod2.reps.imp)


### 
mod2.imp <- lm(rhousevote.x ~ std_survey_ideology+std_survey_ideology_sq+ std_leg_ideology_dems+std_leg_ideology_sq_dems +std_leg_ideology_reps+std_leg_ideology_sq_reps +std_survey_ideology*std_leg_ideology_dems+
std_survey_ideology*std_leg_ideology_reps + factor(pid3) ,data=data2[data2$challengers.total.receipts > 500000,], weight=weight.x)
summary(mod2.comp2)


stargazer(mod2, mod2.comp, mod2.comp2)
